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tiffed in accounting for the spontaneous emission and absorption associated 
' with the coupling to the coherent driving field via an environmental reservoir. 
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driving force - the quantum analog of classical absorption chillers. The depen- 
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dence of the maximum attainable cooling rate on temperature, at ultra-low 
temperatures, is determined and shown to respect the recently-established 
fundamental bound based on the second and third laws of thermodynamics. 
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I. INTRODUCTION 



Recently, a fundamental model-independent bound was derived for the temperature de- 
pendence of the fastest rate Q c at which any substance can be cooled toward absolute zero 
via energetic exchange with a cold reservoir at temperature T c M: 

Qc oc T" , a > 1 . (1.1) 

In the same study, a 3-level quantum-mechanical cooling cycle was postulated and shown to 
achieve this bound. The notion of analyzing molecular cooling with quantum thermodynamic 
cycles had been introduced earlier 0. A system (working fluid) rejects heat to a hot bath, 
removes heat from a cold bath, and is driven by a coherent driving field: an idealized model 
for the laser cooling of magnetically-confined gases at ultra-low temperatures, as well as the 
laser cooling of solids and dyes H|J. 

An important ingredient absent from previous models is the spontaneous emission and 
absorption related to the coupling between the system and the driving field. Superficially, 
this extra coupling provides a dissipative path and hence lowers achievable cooling rate. 
However, we will show that this coupling also introduces the possibility of driving the cool- 
ing cycle thermally, and broadens the conditions under which cooling can be generated. It 
constitutes the quantum-mechanical analog of classical absorption chillers. The quantum 
dynamics and thermodynamic properties of this improved 3-level model will be derived, in 
particular: (a) the identification of thermodynamic with quantum-mechanical variables; (b) 
the conditions required to produce cooling; (c) the dependence of cooling rate on the quan- 
tum control variables; (d) the efficiency of the cooling cycle; and (e) rigorous confirmation 
that the fundamental bound of Eq (|1 . 1| ) is respected as the absolute zero is approached. 

II. THE 3-LEVEL MODEL 

The model is portrayed schematically in Fig [I]: a 3-level system coupled to 3 infinite baths 
(reservoirs) plus an external driving field of coherent radiation. The allowed transitions are: 
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(1) between levels 1 and 3, with a hot bath at temperature Th (heat rejection); (2) between 
levels 1 and 2, with a cold bath at temperature T c (heat removal, i.e., cooling); and (3) 
between levels 2 and 3, simultaneously with the driving field and an environmental bath at 
temperature T e (dissipative heat rejection when T e < Th, and heat input when T e > Th). The 
environmental bath is treated as either: (a) independent of the hot and cold baths; or, when 
interference among the transitions is negligible [f| , (b) representing an additional transition 
to the hot or cold bath. Two special cases are: (1) T e = Th, signifying spontaneous emission 
to the hot bath, i.e., coupling the transition to the hot bath; and (2) T e = T c , describing 
non-radiative (e.g., phonon) decay, i.e., coupling to the cold bath. 
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FIG. 1. Schematic of the 3- level quantum cooling cycle. 
The Hamiltonian of the extended system (system, baths and driving field) is expressed 

as 

H = H s + H sf (t) + H sh + H sc + H se + H h + H c + H e (2.1) 

where H s is the Hamiltonian of the 3-level system (the working fluid of the thermodynamic 
cycle); H s f(t) describes the time- dependent coupling between the system and the driving 
field (coherent radiation); H s h, H sc , H se account for the coupling between the system and 
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its respective baths; and Hh, H c , H e are the respective Hamiltonians of the hot (h), cold (c) 
and environmental (e) baths. 

The system's Hamiltonian can be cast as 

H s = A 3 i P33 + A 21 P 22 (2.2) 

where A y - = Ei — Ej is the transition energy between levels i and j (with E\ chosen as zero), 
and Pij = \i)(J\ for i = j are the projection operators over the states i = 1, 2, 3. (Our units 
are chosen such that h — 1 and fee = 1.) 

Our expression for H s f(t) is based on the semi-classical rotating wave approximation ||: 

H sf {t) = e(P 32 e-^ + P 23 e^) (2.3) 

for a field frequency u and a coupling strength e which depends on the amplitude of the 
driving field as well as the dipole moment of the transition, the latter being assumed to be 
independent of A32. With linear coupling between each transition in the 3-level system and 
its corresponding bath, the other coupling Hamiltonians can be written as: 

H sh = P 31 r h + P 13 r[ , 

H sc = P 21 T C + P 12 Tl , 

H se = P 32 f e + P 23 ft, (2.4) 

where P, L j now represents the raising {i > j) or lowering (i < j) operator for the transition 
% <-> j, and f b denotes a bath operator (b = h,c,e). Only the general properties of the bath 
correlation functions are needed to obtain the reduced description of the system's dynamics. 
The more detailed information on the bath Hamiltonians and the operators Tb will be derived 
afterwards when we examine the explicit forms of the correlation functions. 

The standard tools for quantum open systems can be used to obtain a reduced description 
of system dynamics, i.e., the time evolution for the system's density matrix p s (t), in terms 
of the operators of the 3-level system. p s (t) = Trs{/3(£)} where p is the density matrix of the 
extended system and Tr^ denotes the trace over all bath variables. We invoke the following 
appr oximat ions : 
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• Weak coupling between the 3-level system and each bath. 

• Factorization of the initial density matrix of the extended system, 

p(t = 0) =p a (t = 0)® Ph ® Pc ® Pe (2.5) 
with the initial density matrix of each bath assumed to be in thermal equilibrium, 

g— H b /T b 

p b = — — — (2.6) 
Tr 6 {p 6 j 

• The Markovian assumption that the bath correlations 

C ht , = Tr b {T b f b (-t) p b } = Tr b {T b e^* f' b e* 6 ** p b } (2.7) 

b 

decay rapidly relative to the time scale of evolution of the 3-level system. 

• Weak coupling with the driving field, i.e., e being of the same order as the coupling 
with the baths. This approximation permits us to treat the effect of the field on the 
baths and the system separately 

These approximations result in the density matrix of the extended system factorizing at 
all times. One consequence is that the interfaces between the system and its baths become 
isothermal partitions which allow energy transfer but do not destroy the system's integrity, 
i.e., no quantum entanglement is created between the system and its baths 0. This point 
is essential to the validity of the thermodynamic cycle approach, i.e., to distinct energy 
flows between a system and its reservoirs such that the intactness of the system is not 
compromised. 

The analyses that follow build upon the derivation of the system's reduced dynamics 
as detailed in Appendix [A]. We adopt the Heisenberg representation and express the time 
evolution of an arbitrary operator X as 

BX 

X = i[H s + H sf (t),X] + (— ) + £ C b (X) (2.8) 

b =h,c,e 
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where the super-operator C b describes the effect of the baths on the dynamics of the 3-level 
system and possesses a structure (the Lindblad form) that assures the (complete) positivity 
of the reduced dynamics of Eq ( |2.8|) : 

C b (X) = ^{[P b ,X]Pi - P b [P b \X}} + ^{[Pt,X]P b - P^[P b ,X]} (2.9) 

with P b = P 31 , P 2 i, P 32 for b = h,c,e, repectively. The correlation coefficients and X b , 
defined in Appendix A, depend on the transition energies, the bath temperatures and the 
coupling strength. They appear as coefficients in the equations for P 33 , P 22 , Ai, P32 and 
P23 that follow from Eq ( |2.8| ) and de-couple from the others: 

P 33 = -i e e-^P 32 + i e e^P 23 - (X h + A e )P 33 + A e P 22 + X h P xl 
P 22 =it e-^P 32 - 1 e e^P 23 + A e P 33 - (A c + A e )P 22 + A c P n 
P\x = XhPw + A C P 22 H — (Xh + A c )Pn 

P 32 = 2A 32 P32 - i e e^P 33 + % e e^P 22 - ^{X h + A c + A e + A e }P 32 

P 23 = P 3 f 2 . (2.10) 

The terms P^ are related to the diagonal elements of the reduced density matrix when 
i — j (Tr s {Pjjp s } = pa), and to the off-diagonal elements when i ^ j (Tr s {Py-p s } = p^). 
Hence Eq (|2.10|) represents rate equations for the density matrix elements. The diagonal 
and off-diagonal elements are coupled by the driving field. The bath is responsible for the 
self-couplings of both diagonal and off-diagonal elements. X b and X b are then transition 
probabilities, per unit time, between energy levels. For example, is the transition prob- 
ability of a decay from level 3 to 1 in which energy is rejected to the hot bath, and X^ is 
the probability for an excitation from level 1 to 3 in which energy is absorbed from the hot 
bath. 

III. IDENTIFICATION OF THE THERMODYNAMIC VARIABLES 

The identity of the thermodynamic energy flows in the cooling cycle follows from energy 
conservation and averaging The Hamiltonian Ht = H s + H s f(f) is introduced into the 
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evolution equation (|2.8|) , is multiplied by the initial density matrix, and is traced over the 
variables of the 3-level system, to yield 

(4) = <^^> + {£h(Ht)) + (£c(Ht)} + (C e (H T )) . (3.1) 

The energy flow associated with the driving field (the first term on the RHS of Eq ( |3.1| )) 
is the power input (cycle-averaged work) to the cycle, W. The remaining 3 terms are the 
respective heat flows between the system and its baths, Qh, Qc, Qe- At steady state, (Ht) 
is constant and independent of the system's initial state, so (Ht) = and Eq ( |3.1|) can be 
expressed as the First Law of thermodynamics for the thermodynamic cycle: 

W + Q h + Q c + Q e = (3.2) 

with energy flows into the system defined as positive. 

Now we can derive the relation between the cycle's thermodynamic variables and the 
quantum-mechanical parameters. At steady state (denoted by the superscript ss), the gen- 
eral P operators become 

p ij -> e ia ^ (3.3) 

where P^ = 0, an = 0, a^2 = (jJ and the remaining values of P?? follow from introducing 
( p.3|) into Eq ( p.lOj) and solving the resultant set of linear time-independent equations. The 
solutions for P^ turn out to be proportional to the identity operator 1, P?? = Pijl, where 
is real when i = j and complex when i ^ j. From the relation between the P operators and 
the elements of the density matrix, the steady-state expectation values of the populations 
(i.e., the diagonal elements) are constants, (Pa) = (Pff) = pa, and the expectation values 
of the off-diagonal elements P32 and P23 oscillate in time, (Pij) = (Pff e iaijt ) = p^ e taijt . 
Explicit formulae for are provided in Appendix [B] for the condition of resonance that is 
assumed in the analysis that follows 

00 = A 32 . (3.4) 
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The cycle-average thermodynamic variables of the quantum refrigerator can now be 
expressed as: 

W = 2 e A 32 Imag[p 32 ] = A 32 e 2 {A(n c - n h ) + B} 

Q h = A 3 i (-\hP33 + XhPu) = -A 31 {e 2 A(n c - n h ) + C(n c n e - n h )} 

Q c = A 2 i (-X c p 22 + XcPu) = A 21 {e 2 A(n c - n h ) + C(n c n e - n h )} 

Q e = A 32 (-\eP33 + KP22) = -A 32 {e 2 B - C(n c n e - n h )} (3.5) 

where Eq ( |A6|) of Appendix ^ has been used, and n b denotes the equilibrium populations: 
rih = e~ Am ^ Th , n c = e~ A21//Tc , n e = e _A,32 / Te . Detailed expressions for the positive coeffi- 
cients A,B and C are presented in Appendix B, and depend on the bath temperatures, the 
transition energies and the coupling strength. Each energy flow constitutes a competition 
between two processes with rates A& XbPkk and A& XbPu (k > I), with the former describing 
heat rejection of magnitude A& to the bath in the decay k — > I, and the latter representing 
heat removal of the same magnitude from the bath in the excitation I — > k. Also, as will 
be elucidated in Section 0, the work input rate W need not be the only driving force for 
the cooling cycle. When T e > Th, the incoherent thermal flow Q e can become positive and 
hence can also contribute to the cooling rate. 
The entropy production rate for the cycle is 

J-h J-c J-e 

= -4e 2 ^(n c - n h )(^ - ^L) + Ae ^ + X h ) A^ (1 _ + 

XhX c X e . \/^3i A21 A 32 . , . 
jj— (nh-n c n e )(— — ) (3.6) 

where D is a positive function derived in Appendix EL With the relations noted above for 



the equilibrium populations rib, it is straightforward to prove that S in Eq (|3.6| ) must be 
non-negative, in accordance with the Second Law. 
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IV. COUPLING WITH THE BATHS 



The results for the cooling performance of the quantum refrigeration cycle depend on the 
particulars of the coupling between the system and its baths. Hence specific models must 
be invoked. The primary variable of interest is the cooling rate (the interaction with the 
cold bath). Since in any event the results that will now be derived are independent of the 
nature of the hot and environmental reservoirs, we will treat these two as white baths, i.e., 
baths with a constant density of energy modes (at least in the frequency range of interest), 
for which 

X b = A b , \ b = A b e~ A »/ T > (4.1) 

where A& is the strength of the coupling. 

Earlier quantum refrigeration models also treated the cold bath as having a constant 
mode density |l]J^]. Here a more elaborate and realistic model is introduced for the cold 



bath: an assembly of harmonic oscillators - a viable model in the weak coupling limit ||T0 
The bath's Hamiltonian can then be expressed as 

H c = ]T A Ci al 4 a C4 (4.2) 

i 

where a C j and a ci are the annihilation and creation operators, with the index % spanning the 
bath oscillator energies A C j. The operator of Eq ( |2.4| ), that couples the system and bath is 
then 

r c = J2 9cta C i + g* ci a ] ci (4.3) 

i 

where g ci and g* ci denote the coupling constants. 

The influence of an harmonic bath on system dynamics is embedded in the spectral 
strength function |10|JTT[| J C (A) = J2i \9ti\ 2 $(A — A ci ) which, for low energies, is well approx- 



imated by a simple power law dependence 

J C (A) « A c A 5 ' (4.4) 



where s c must be positive [I2[. The coefficients A then follow as 



A c ( A c , T c 



1 _ e -A c /T c 



A C (A C , T c 



3 A C /T C _ J 



(4.5) 



V. COOLING PERFORMANCE 
A. Objectives 

Several aspects of the thermodynamic performance of the quantum cooling cycle will be 
addressed in this section: (1) how strong a coupling with the driving field is required to 
achieve cooling (i.e., to insure Q c > ); (2) how does cooling rate vary with the natural 
control variable A 2 i and under what conditions is cooling rate maximized; (3) how does 
maximum cooling rate depend on T c in the limit of the absolute zero and how does this 
dependence compare to the fundamental bound established from the Second and Third 
Laws; and (4) what is the efficiency of this cooling cycle and can it be cast in the same form 
as classical chiller analyses. In this section, only T e values in the range T c < T e < are 
analyzed, so the environmental bath serves a solely dissipative role. The following section 
will focus on the special thermodynamic consequences when T e > T^. 

B. Cooling window and electronic analog 

Because of the competition between the coupling to the driving field and dissipative 
losses, only certain combinations of system parameters produce cooling. The thermoelectric 
chiller offers a familiar example of the existence of a refrigeration window, where the cooling 
effect exists only for: (a) a sufficiently high thermopower (Seebeck coefficient); and (b) a 
particular voltage window. There are analogous limits in the 3-level quantum chiller. 

In fact, the 3-level model can be viewed as an electronic device, with the analog of the 
voltage V being A 32 (in units of electronic charge). The maximum voltage is V maiX = A 31 , 
and A 2 i = V max — V. The power input is the product of voltage and current, so that 
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Eq (|3.5|) provides the analog of the electrical current as a complicated implicit function of 
the system parameters, as well the relation between current and voltage. In the absence 
of the environmental bath, current is an exponentially increasing function of voltage |TJ 
which, not coincidentally, is the same as for ideal diodes. With the added dissipation to the 
environmental bath, the current remains a strongly increasing function of voltage, albeit not 
strictly exponential, similar to non-ideal diodes. 

The minimum coupling strength with the driving field e m i n to produce a cooling effect 
follows from Eq ( ^5|) : 

e m in = \J C(n h - n e )/{A(l - n h )} (5.1) 
and is plotted clS 8u function of T e in Fig 0. 
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FIG. 2. Minimum coupling needed for cooling, as a function of the environmental bath temper- 
ature, at 3 values of coupling coefficient with the environmental bath. Other system parameters 
are fixed at: T c = 0.1, = 0.2, A31 = 1, = A c = 0.001, and s c = 1. The regime of T e > T^, 



where no coupling is required to produce cooling, is addressed in Section VI 
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C. Maximum cooling rate and the fundamental bound 

Our "control knob" for varying cooling rate is the transition energy A 21 . The cooling rate 
vanishes at two values of A 21 that delimit the cooling window: A 21 = and A 2 i = AJJf"*, 
where A™ x depends on the principal system parameters, most notably T c and A e (Eq (|3.5|) ), 
as illustrated in Fig |3|. Sample curves of cooling rate as a function of A 21 are plotted in Fig 




FIG. 3. Refrigeration window as a function of T c for 3 values of the coupling coefficient with 
the environmental bath. Other system parameters are fixed at: = T e = 0.2, A31 = 1, 
A h = A c = e = 0.001 and s c = 1. 
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FIG. 4. Cooling rate as a function of A21 for 3 values of cold bath temperature. Other system pa- 
rameters are fixed at: Ty l = T e = 0.03, A31 = 1, = A c = A e = e = 0.001 and s c = 1. The broken 
curve and solid circles indicate the locus of maximum cooling rate. The solid curves are nominally 
exact numerical results, while the solid squares illustrate the accuracy of the approximation of Eq 



Solving for the maximum cooling rate numerically in the limit T c — > reveals that this 
3- level quantum chiller respects the fundamental bound of Eq ( |1 . 1| ) . A more satisfying 
analytic derivation is possible once we establish the accuracy of an approximate closed-form 
expression for the cooling rate in the limit of the absolute zero. In the limit of vanishing T c : 
(a) A^ ax vanishes; (b) A 32 w A 31 ; and (c) the coefficients c w , c 20 and c^q of Appendix [B| 
become independent of both T c and A 21 . Then, to an excellent approximation, the cooling 
rate reduces to 



^21 



A, 



(e -A 21 /T c Ciq _ C2q) 



(5.2) 



1 _ e -A 21 /T c CdQ 

The accuracy of this approximation is depicted in Fig |j. Although maximizing the cooling 
rate of Eq ( |5.2| ) yields a transcendental equation that can only be solved numerically, the 
solution for that transition energy A 21 is of the form A 21 oc T c , independent of the coupling 
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with the environmental bath or the value of the positive exponent s c in Eq ( |4.4| ). It then 
follows from Eq ( |5.2| ) that 

gmax ^ T (, c +1) ^ > (53) 

as required by the Second and Third laws ffl. 

D. Efficiency 




cooling rate 

FIG. 5. Characteristic chiller plot, COP against cooling rate, with and without coupling to the 
environmental bath. The control variable is A21. System parameters are fixed at: T c = 0.1, 
T h = T e = 0.2, A31 = 1, A h = A c = e = 0.001 and s c = 1. 

Cooling cycle efficiency is usually defined by Coefficient of Performance (COP), the ratio 
of cooling rate to input power, which in this instance is COP = Q c /W. Cooling cycles are 
conveniently characterized by a plot of COP against cooling rate ]nj , as in Fig |5]. Even in 
the absence of the (parasitic) environmental bath, the 34evel system possesses an energy 
leak that militates against efficient operation as A21 — > Q, which appears as the lower 
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branch of the curves in Fig |5], with COP vanishing as cooling rate is lowered. At the other 
end of the refrigeration window ( A 2 i — > A^ ax ) , the existence of the extra dissipation to the 
environmental bath makes a quantum difference. In the absence of this extra energy-leak 
pathway, there are no irreversibilities that mitigate against efficiency operation, so both 
cooling rate and power input vanish at the same rate such that the COP approaches its fun- 
damental reversible (Carnot) value of T c j (Th — T c ) |I|Jl3|. The coupling to the environmental 
bath introduces a loss mechanism that mitigates against this nominally slow operation, so 
COP also vanishes in this limit (the upper branch of the curve in Fig ||for A e > 0). Fig |6| 
offers an alternative view of the differences in dissipation: a plot of entropy production rate 
against A 21 with and without coupling to the environmental bath. 




FIG. 6. Entropy production rate as a function of A21, with and without coupling to the en- 
vironmental bath. The refrigeration window is broader in the absence of this coupling (the two 
curves terminate at different values of A^*). Arrows indicate the points of maximum cooling rate. 
System parameters are fixed at: T c = 0.1, = T e = 0.2, A31 = 1, = A c = e = 0.001 and 
s c = 1- 
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VI. QUANTUM ABSORPTION CHILLER 



Classical cooling cycles can be driven by thermal sources that are hotter than the heat 
rejection reservoir, commonly called absorption chillers (in contrast to more common work- 
driven mechanical chillers) ||13|| . A simple quantum analog is our 3- level cooling cycle with 
T e > Tfr, where the environmental bath is analogous to what is called the generator in clas- 
sical absorption cycles. Consider de-coupling the 3-level system from the coherent radiation 
field and driving it with incoherent radiation (e.g., a flash lamp) at input thermal power 
Q e . The environmental bath becomes a heat source rather than a dissipative sink. As in 
classical absorption chillers, the COP is defined as Q c / Q e . While the refrigeration window 
and characteristic chiller curve will now be derived, we note that the basic result for the 
dependence of maximum cooling rate on T c in the limit T c — > ( |5.3|) remains unaltered. 

The refrigeration window follows from Eq Q3.5|) : 



A max T c (T e Th) , , 

A 21 = TfTTTf TfT A 31 (6.1) 

J-h \J-e — J-c) 

which, as for the dissipative environmental bath at T e < Th, vanishes in the limit T c — > 0. 

The characteristic chiller curve is graphed in Fig [7[ With T e > Th, there is no irreversibil- 
ity that undermines efficient slow operation ( A 21 — > A^ ax ) . Hence the ordinate intercept of 
each curve can approach the reversible Carnot limit for absorption cycles of (TB| 

j i_ 

COPcarnot = ~f f~ ■ (6.2) 
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FIG. 7. COP against cooling rate for the 3-level quantum absorption cycle, at 3 values of T e . 
System parameters are fixed at: T c = 0.1, Th = 0.2, A31 = 1, A^ = A c = A e = e = 0.001 and 
s r = 1. 



Sections [V| and |VI] addressed the extreme cases of the power input to the cooling cycle 
being either pure coherent radiation or purely incoherent and thermal, respectively. A 
combination of the two (when T e > Th) can be calculated directly from the results derived 
above. 



VII. SUMMARY 

In approaching the absolute zero, any cooling cycle must be dominated by quantum 
dynamics. At the same time, one should be able to probe its behavior with fundamental 
chiller thermodynamics. The 3-level quantum model proposed and explored here provides a 
simple case study. It steps beyond earlier quantum refrigeration models in accounting fully 
for spontaneous emission (and spontaneous absorption), and including an environmental 
bath that either mimics actual rate-dependent dissipative mechanisms in work-driven chillers 
when T e <Th, or establishes the quantum analog of an absorption (purely thermally-driven) 
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chiller when T e > 7\. All the irreversibilities modeled here and in previous studies |],|2| 
are effectively heat leaks. The quantum analog of internal friction stems from the non- 
commutation of the field and system Hamiltonians. In the limit of week fields (e — > 0), this 
internal dissipation grows negligible. Incorporating the quantum version of friction remains 
a challenge for a future study. 

Subject to reasonable approximations that retain the integrity of the system and its 
reservoirs, expressions can be derived for: (a) the conditions under which a cooling effect 
can be generated, (b) chiller efficiency as a function of cooling rate, and (c) the temperature 
dependence of the fastest rate at which molecular systems can be cooled at ultra-low tem- 
peratures. The relation of this maximum cooling rate to the energy mode density has been 
established, and the basic result has been shown to respect the fundamental bound that was 
recently established for maximum cooling rate from the Second and Third Laws. 
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APPENDIX A: QUANTUM REDUCED DYNAMICS 

To obtain the evolution equation for the 3-level system, we invoke standard tools for 
quantum open systems: the Redfield approach |I^JT5|| and the secular approximation |IT| , [16| . 



Assuming (a) weak coupling of the system with its bath, (b) a weak driving field, (c) rapidly 
decaying bath correlation functions, and (d) an initial density matrix in the form of Eq Q2.5|) , 
we can express the evolution equation as: 

'p s (t) = -i[H s + H sf {t),p s (t)) + 

{-GWA31) [P31, Pi3Ps(t)) - G ftfh (-A 31 ) [As, + {H.C.} + 

™ h h tl 

{-G t ft (A 21 ) [Ai, Pi2Ps(t)] - G ftf ,(-A 21 ) [A2, AiP.(t)]} + {H.C.} + 



{-G fef t(A 32 ) [P32, P23Ps(t)} - % fe (-A 32 ) [As, P32Ps(t)}} + {H.C.} 

(Al) 

where {H.C.} denotes the Hermitian conjugate of the expression in brackets that immedi- 
ately precedes it. Eq ( |AI] ) is valid independent of the nature of the baths as long as the 
coupling is weak and the bath correlation functions decay quickly. 
The influence of the baths is included in the coefficients G, 

/•CO 

G ht ,(A b ) = j Q dte UA »C tbt , (A2) 

with A b = A 3 i, A 2 i, A 32 for b = h,c,e, respectively, and Cf b f, being given by (|2.7| ). When 
f' b = f I, G can be decomposed as 

Cppt(A 6 ) 

G tbtl (A b ) = —Lt + iS hn (A b ) (A3) 

with C, S being real and C > [gJT^JT^] . With Eq (|A3|), we split each term of Eq (|Al|) into 



two parts related to C and S. The C-related term is 

C tht i(A b ) A . C t , t (-A b ) 

- — + \p.m>p}]} - — {tn T , Ap s (t)] + [p fl (t)pt, p b]} 

(A4) 

with P b = P31, Ai, -P32 for b = h,c,e respectively. It is readily confirmed that the terms of 
have the standard Lindblad form that insures complete positivity of the dynamics. 
For simplicity of notation in Sections |Il|-[lV| and Appendix [B], we introduce the notation 



\ b = \ b (A b ,T b ) = C tht ,(A b ) , X b = A 6 (A fe ,T b ) = C m (-A b ) (A5) 

" b b v 

with X b and \ b related by 

X b = e -^/n Xb . (A6) 

The S'-related term is 

iS t bP (A b )[P b P b \p s (t)} +iS t i t (-A b )[PiP b ,p s (t)} . (A7) 

b h 
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PbPb and PI Pb represent small corrections to the system's energy levels, which in the weak 
coupling limit also turn out to be negligible. 

Hence, from Eq (|Al|), the time evolution equation for any operator (in the Heisenberg 
representation) can be obtained, and the germane results are provided in Section ||. 

APPENDIX B: COEFFICIENTS AT STEADY STATE 

The coefficients Pij are the solution of a 4 x 4 system of linear equations the coefficients 
of which depend on: (a) the system energy structure, (b) the coupling with the field e, and 
(c) the coupling with the baths and 

Pn = (cio + en A c + C12 X 2 C )/D 

P22 = (c 20 + c 2 i A c + c 2 i A c + c 22 X C X C )/D 

P33 = (c 3 o + c 3 i A c + c 3 i A c + c 32 A 2 + c 32 XAc)/D 

P32 = i (C c0 + C c i A c + Cd \ C )/D , P23 = p* 32 (Bl) 

where the denominator D is 

D = c d0 + c d i A c + c d i A c + c d2 X 2 C + c d2 A C A C • (B2) 

and the c coefficients are given by 

c d0 = (X h + X e + X e ){{X h + X h )X e + X e X h } + Ae 2 {X h + 2X h ) 

c d i = (X h + A e + A e ) (A ft + X h + A e ) + {X h + X h )X e + X e X h + 4e 2 

c d i = (X h + Ac + Ac) 2 + 8e 2 

Cd2 = Xh + Xh + A e 
Cd2 = Xh + A e + A e 

cio = (X h + A e + Ae)AftAe + 4e 2 A ft 

en = (X h + Ac + Ac) (Aft + A e ) + Aft A e + 4e 2 

C12 = (Aft + A e ) 
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C20 


= (Xh + A e + A e )A e A h + Ae 2 X h 


C 2 1 


= A e Afe 


P21 


= (A/, + A e + A e )(A h , + A e ) + 4e^ 


C22 


= (Aft + A e ) 


c 30 


= (A ft + A e + A e )A e A/j + 4e A/> 


C 3 1 


= (Aft + A e + A e )Ah + A e Aft 


C 3 1 


= (A h + A e + A e )A e + 4e 2 


C32 


= A ft 


C32 


= A e 


C C 


= 2eAft(A e - A e ) 


C C 1 


= -2eAft 


C C 1 


= 2e(Aft + A e -A e ). 



The coefficients A, B and C of Eq flO]) can be expressed as 

A = Ae 2 \ h \ c /D 

B = A{\ h + \ c ){\ e -\ e )/D 

C = AftA c A e (Aft + A c + A e + X e )/D . (B4) 
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